{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "0NQT88nVbsOv"
   },
   "source": [
    "# 1. Описание набора данных и признаков"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "kCw3EejCbsOy"
   },
   "source": [
    "Данные взяты с портала http://stat.gibdd.ru/\n",
    "\n",
    "Для получения данных применен парсинг ответов REST API, который используется для фронт-энда вышеуказанного сайта. Информация собрана со всех регионов за последние три года. Код парсера можно найти в [файле parser.ipynb](https://drive.google.com/file/d/1mjqL7aKz36l6qtHc8igVNNdmWCddE9Bt/view?usp=sharing).\n",
    "В нем мы опрашиваем бэкэнд сайта, по всем существующим регионам и районам (список ОКАТО) за каждый месяц, складывая ответ в отдельные файлы. Я выполнил сохранение исходного датасета в базу данных Mongo с последующим переводом в сжатый csv.\n",
    "\n",
    "Скачать исходный датасет можно тут: [Ссылка](https://drive.google.com/open?id=1KNTJz-peQ4vSzSgTcuXaXu6R4Fygvwu8)\n",
    "\n",
    "Ссылка для скачивания данных о населении регионов: [Ссылка](https://drive.google.com/file/d/1D-yoaTU4Zjx3h-iouMKyZklGBKGlw35X/view?usp=sharing)\n",
    "\n",
    "Задача - понять причины, наиболее характерные для опасных ДТП, приведшие к гибели участников, их относительную важность. Провести анализ, обучить подходящую модель и получить веса признаков, которые и покажут важность факторов.\n",
    "Данную модель и вес признаков можно использовать для себя, при вождении, для понимания, где лучше ездить помедленнее и поаккуратнее. Или автопилотах, с той же целью, внимательнее рядом с бетономешалками.\n",
    "\n",
    "На сайте представлены только аварии, где был как минимум один раненый или один погибший, то есть записей без пострадавших нет вообще, поэтому **как целевой признак выбран признак наличия погибших**.\n",
    "\n",
    "\n",
    "Каждая строка датасета (всего строк 501367) представляет из себя JSON объект с параметрами ДТП примерно такого вида:"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "LLe7PucabsOz"
   },
   "source": [
    "```python\n",
    "{'COORD_L': '38.052969',  # координаты\n",
    " 'COORD_W': '55.631653',\n",
    " 'DTP_V': 'Наезд на пешехода',  # вид дтп\n",
    " 'District': 'Раменский р-н',  # район\n",
    " 'K_TS': 1,  # количество транспортных средств, учавстовавших в дтп\n",
    " 'K_UCH': 2,  # количество людей-участников\n",
    " 'KartId': 198310083,  # ИД объекта в базе ГИБДД\n",
    " 'OBJ_DTP': ['Остановка общественного транспорта',  # объекты интереса рядом\n",
    "  'Нерегулируемый пешеходный переход',\n",
    "  'Жилые дома индивидуальной застройки',\n",
    "  'Гаражные постройки (гаражный кооператив, товарищество либо иное место концентрированного размещения гаражей)'],\n",
    " 'POG': 0,  # число погибших, целевой признак!\n",
    " 'RAN': 1,  # число раненых\n",
    " 'Time': '08:20',  # время ДТП\n",
    " '_id': 198310083,  # id для Mongo, совпадает с KartId\n",
    " 'change_org_motion': 'Режим движения не изменялся',  # последствия ДТП, частичное или полное перекрытие\n",
    " 'date': '31.08.2016',  # дата\n",
    " 'date_name': 'авг 2016',  # она же текстово\n",
    " 'dor': 'МОСКВА-ЖУКОВСКИЙ-А/ПБЫКОВО',  # дорога\n",
    " 'dor_k': '6',  # категория дороги, качество\n",
    " 'dor_z': 'Региональная или межмуниципальная (дорога регионального или межмуниципального значения)',  # тип\n",
    " 'factor': ['Сведения отсутствуют'],  # дополнительные факторы\n",
    " 'filename': '46_46248_2016_08.json',  # имя файла парсера\n",
    " 'house': '',  # дом\n",
    " 'k_ul': '',  # район\n",
    " 'km': '2',  # положение на дороге км и м дальше\n",
    " 'm': '100',  # \n",
    " 'n_p': 'дп Удельная',  # населенный пункт\n",
    " 'ndu': ['Отсутствие, плохая различимость горизонтальной разметки проезжей части',  # дорожные условия\n",
    "  'Отсутствие дорожных знаков в необходимых местах',\n",
    "  'Неправильное применение, плохая видимость дорожных знаков'],\n",
    " 'osv': 'Светлое время суток',  # освещение\n",
    " 'reg_name': 'Московская область, Раменский район',  # регион, район\n",
    " 'rowNum': 1,  # номер строки для фронт-энда\n",
    " 's_dtp': '820',  # тип схемы ДТП, по нему можно смотреть картинки тут http://stat.gibdd.ru/img/SDTP/820.png , где последняя цифра в имени файла - этот тип\n",
    " 's_pch': 'Мокрое',  # состояние дороги\n",
    " 's_pog': ['Пасмурно', 'Дождь'],  # погода\n",
    " 'sdor': ['Перегон (нет объектов на месте ДТП)', 'Иное место'],  # объекты непосредственно на месте ДТП (в отличие от OBJ_DTP, где указаны объекты поблизости)\n",
    " 'street': '',  # улица\n",
    " 'ts_info': [{'color': 'Серый',  # цвет машины\n",
    "   'f_sob': 'Частная собственность',  # вид собственности\n",
    "   'g_v': '1995',  # год выпуска\n",
    "   'm_pov': '',  # повреждения\n",
    "   'm_ts': 'Прочие модели КАМАЗ',  # модель\n",
    "   'marka_ts': 'КАМАЗ',  # марка\n",
    "   'n_ts': '1',  # номер тс\n",
    "   'o_pf': 'Физические лица',  #\n",
    "   'r_rul': 'Задний',  # тип привод\n",
    "   't_n': 'Технические неисправности отсутствуют',  # тех неисправности\n",
    "   't_ts': 'Прочие',  # класс автомобиля, \n",
    "   'ts_s': 'до 1 сут.',  # ?\n",
    "   'ts_uch': [{'ALCO': '',  # алкоголь в крови, число промилле\n",
    "     'INJURED_CARD_ID': '',  # карточка медицинская\n",
    "     'K_UCH': 'Водитель',  # категория, водитель или пассажир\n",
    "     'NPDD': ['Нет нарушений'],  # нарушение ПДД данным участником\n",
    "     'N_UCH': '2',  # номер для фронтэнда\n",
    "     'POL': 'Мужской',  # пол\n",
    "     'SAFETY_BELT': 'Да',  # ремень бехопасности\n",
    "     'SOP_NPDD': ['Оставление места ДТП'],  # сопутсвующие нарушения ПДД\n",
    "     'S_SEAT_GROUP': '',  # ряд сидений ?\n",
    "     'S_SM': 'Скрылся, впоследствии разыскан (установлен)',  # статус\n",
    "     'S_T': 'Не пострадал',  # травмы\n",
    "     'V_ST': '10'}]}],  # водительский стаж\n",
    " 'uchInfo': [{'ALCO': '',\n",
    "   'K_UCH': 'Пешеход',\n",
    "   'NPDD': ['Переход через проезжую часть вне пешеходного перехода в зоне его видимости либо при наличии в непосредственной близости подземного (надземного) пешеходного перехода'],\n",
    "   'N_UCH': '1',\n",
    "   'POL': 'Женский',\n",
    "   'SOP_NPDD': ['Нет нарушений'],\n",
    "   'S_SM': 'Нет (не скрывался)',\n",
    "   'S_T': 'Раненый, находящийся (находившийся)  на амбулаторном лечении, либо которому по характеру полученных травм обозначена необходимость амбулаторного лечения (вне зависимости от его фактического прохождения)',\n",
    "   'V_ST': ''}]}\n",
    "```"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "kQ_E97WNbsO1"
   },
   "outputs": [],
   "source": [
    "import numpy as np\n",
    "import pandas as pd\n",
    "pd.options.display.max_columns = 12\n",
    "# Disable warnings in Anaconda\n",
    "import warnings\n",
    "warnings.simplefilter('ignore')\n",
    "# We will display plots right inside Jupyter Notebook\n",
    "%matplotlib inline\n",
    "import matplotlib.pyplot as plt\n",
    "# We will use the Seaborn library\n",
    "import seaborn as sns\n",
    "sns.set()\n",
    "# Graphics in SVG format are more sharp and legible\n",
    "# %config InlineBackend.figure_format = 'svg'\n",
    "\n",
    "\n",
    "import os\n",
    "import json\n",
    "import pprint\n",
    "from tqdm import tqdm\n",
    "from itertools import chain\n",
    "\n",
    "from sklearn.preprocessing import MultiLabelBinarizer, OneHotEncoder, LabelBinarizer\n",
    "from sklearn.linear_model import LogisticRegression,  LogisticRegressionCV, SGDClassifier\n",
    "from sklearn.ensemble import RandomForestClassifier, ExtraTreesClassifier, GradientBoostingClassifier, AdaBoostClassifier\n",
    "from sklearn.model_selection import cross_val_score, GridSearchCV, train_test_split, validation_curve, learning_curve, KFold\n",
    "from sklearn.svm import LinearSVC, SVC\n",
    "from sklearn.neighbors import KNeighborsClassifier\n",
    "from sklearn.metrics import roc_auc_score\n",
    "\n",
    "import eli5\n",
    "\n",
    "\n",
    "RANDOM_STATE = 17"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "HCAZEXZWhPa2"
   },
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "CXD-GyMSbsO7"
   },
   "source": [
    "# 2. Первичный анализ признаков"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "8yHw2-kIcgFi"
   },
   "outputs": [],
   "source": [
    "data = pd.read_csv('gidbb_info.tar.gz')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     },
     "base_uri": "https://localhost:8080/",
     "height": 34
    },
    "colab_type": "code",
    "executionInfo": {
     "elapsed": 544,
     "status": "ok",
     "timestamp": 1524235623538,
     "user": {
      "displayName": "Павел Куликов",
      "photoUrl": "//lh6.googleusercontent.com/-hT0LtsqJaKc/AAAAAAAAAAI/AAAAAAAAD9k/IHmlCLzzyMM/s50-c-k-no/photo.jpg",
      "userId": "106446428922776960365"
     },
     "user_tz": -180
    },
    "id": "yENT8wnabsPC",
    "outputId": "7d83e40a-d3de-4beb-cd6e-11cb4acda819"
   },
   "outputs": [],
   "source": [
    "data.shape"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "gZLmICPabsPM",
    "outputId": "0fc24d48-46a7-487d-dca0-e97b534c9877"
   },
   "outputs": [],
   "source": [
    "data.head()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "EQIH9irZbsPT"
   },
   "outputs": [],
   "source": [
    "# убираем совсем ненужные поля\n",
    "\n",
    "data.drop(['_id', 'km', 'm', 'KartId', 'rowNum'], inplace=True, axis=1)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "ROB8WnhxbsPX",
    "outputId": "47d8aeda-b95e-4f0b-a30d-a0f043972caa"
   },
   "outputs": [],
   "source": [
    "data.describe()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "miFUyWI9bsPd",
    "outputId": "892cdcf8-0842-4319-ba37-4539a2c39403"
   },
   "outputs": [],
   "source": [
    "data.info()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "JlENrtj0bsPl"
   },
   "source": [
    "Видно, что есть пропуски. Колонки n_p, street, house, dor, change_org_motion использовать в дальнейшем не будем, это просто названия городов, дорог и т.п.\n",
    "Есть две интересующие нас колонки с пропусками, посмотрим их подробнее."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "Oghqmm-hbsPm",
    "outputId": "d75dcad1-5b20-4bc6-95c5-65ff3ae17a9e"
   },
   "outputs": [],
   "source": [
    "data.k_ul.value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "cANLGC1RbsPs",
    "outputId": "66edc73e-59ec-44f3-e9ae-4a4b7122186b"
   },
   "outputs": [],
   "source": [
    "data.dor_k.value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "nLVOncvMbsP0"
   },
   "outputs": [],
   "source": [
    "# в первом случае заменяем отсутвующие значения на отдельное значение\n",
    "data.k_ul.fillna('Не указано', inplace=True)\n",
    "# во втором можно бы было заменить на медиану, поскольку признак отражает качество дороги, 1 наилучшее, 7 наихудшее, по сути направление. Но имеющихся данных слишком мало, треть датасета, лучше заменим на что-то новое, всё расматривать будем как категорию\n",
    "data.dor_k.fillna(0, inplace=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "kd-_I3rwbsP9"
   },
   "outputs": [],
   "source": [
    "# раскроем json в полях верхнего уровня\n",
    "for c in ['OBJ_DTP', 'factor', 'ndu', 's_pog', 'sdor']:\n",
    "    data[c] = data[c].apply(json.loads)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "A86P4k1ObsQB"
   },
   "outputs": [],
   "source": [
    "data['ts_info'] = data['ts_info'].apply(json.loads)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "kFD3gZswbsQH"
   },
   "outputs": [],
   "source": [
    "# Дозаполним поля датафрейма из большого словаря ts_info\n",
    "def get_sublists_npdd(x):\n",
    "    res = list()\n",
    "    for rec in x:\n",
    "        for ts in rec['ts_uch']:\n",
    "            for npdd in ts['NPDD']:\n",
    "                res.append(npdd)\n",
    "    return res\n",
    "                    \n",
    "data['t_ts'] = data.ts_info.apply(lambda x: [elem['t_ts'] for elem in x])\n",
    "data['color'] = data.ts_info.apply(lambda x: [elem['color'] for elem in x])\n",
    "data['V_ST'] = data.ts_info.apply(lambda x: [i['V_ST'] for sublist in x for i in sublist['ts_uch']])\n",
    "data['SAFETY_BELT'] = data.ts_info.apply(lambda x: [i['SAFETY_BELT'] for sublist in x for i in sublist['ts_uch']])\n",
    "data['S_T'] = data.ts_info.apply(lambda x: [i['S_T'] for sublist in x for i in sublist['ts_uch']])\n",
    "data['all_NPDD'] = data.ts_info.apply(get_sublists_npdd)\n",
    "data['ALCO'] = data.ts_info.apply(lambda x: [i['ALCO'] if i['ALCO'] else 0 for sublist in x for i in sublist['ts_uch']])\n",
    "data['POL'] = data.ts_info.apply(lambda x: [i['POL'] for sublist in x for i in sublist['ts_uch']])\n",
    "data['m_ts'] = data.ts_info.apply(lambda x: [elem['m_ts'] for elem in x])\n",
    "data['marka_ts'] = data.ts_info.apply(lambda x: [elem['marka_ts'] for elem in x])\n",
    "data['t_n'] = data.ts_info.apply(lambda x: [elem['t_n'] for elem in x])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "uwtCSEnYbsQL",
    "outputId": "6e3fb75d-591d-45c5-cf2e-5308bb6fd44b"
   },
   "outputs": [],
   "source": [
    "data.POG.value_counts().plot(kind='bar', figsize=(12,4))\n",
    "plt.legend()\n",
    "plt.title('Распределение погибших');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "YrVL7uQlbsQU",
    "outputId": "704c9490-4573-4576-aeba-00e9504e27f1"
   },
   "outputs": [],
   "source": [
    "data.RAN.value_counts().plot(kind='bar', figsize=(12,4))\n",
    "plt.legend()\n",
    "plt.title('Распределение раненых');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "aQsFz8wAbsQb",
    "outputId": "c0c7b7f4-afb1-45fb-bcb1-1c45b21a8dae"
   },
   "outputs": [],
   "source": [
    "#  строк, где никто не погиб и не ранен, нет\n",
    "len(data[(data.RAN == 0) & (data.POG == 0)])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     },
     "base_uri": "https://localhost:8080/",
     "height": 500
    },
    "colab_type": "code",
    "executionInfo": {
     "elapsed": 654,
     "status": "ok",
     "timestamp": 1524235642229,
     "user": {
      "displayName": "Павел Куликов",
      "photoUrl": "//lh6.googleusercontent.com/-hT0LtsqJaKc/AAAAAAAAAAI/AAAAAAAAD9k/IHmlCLzzyMM/s50-c-k-no/photo.jpg",
      "userId": "106446428922776960365"
     },
     "user_tz": -180
    },
    "id": "D_EsTZxJbsQl",
    "outputId": "ddb16054-3cc6-4b24-b855-6c717ec32fe4"
   },
   "outputs": [],
   "source": [
    "# создаем целевую переменную\n",
    "data['target'] = data['POG'] > 0\n",
    "\n",
    "data.target.value_counts().plot(kind='bar', label='target')\n",
    "plt.legend()\n",
    "plt.title('Распределение целевого признака');"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data.target.value_counts()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "HC7ABjBQbsQ1",
    "outputId": "b4bf484f-f221-43ed-a30f-80459e3276b9"
   },
   "outputs": [],
   "source": [
    "# большинство признаков - категориальные, поэтому корреляции трудно оценить на этапе первичного анализа\n",
    "sns.heatmap(data.corr())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Естественна корелляция между таргетом и числом погибших, легкая корреляция между погибшими и ранеными, сильная корреляция между количеством машин/участников/раненых, что понятно."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "YedBTjM9bsQ6"
   },
   "source": [
    "# 3. Первичный визуальный анализ данных"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ct = pd.crosstab(data['DTP_V'], data['target'], normalize=True)\n",
    "\n",
    "ct.plot.bar(stacked=True)\n",
    "plt.legend(title='target')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ct = pd.crosstab(data['change_org_motion'], data['target'], normalize=True)\n",
    "\n",
    "ct.plot.bar(stacked=True)\n",
    "plt.legend(title='target')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data['SAFETY_BELT_NO'] = data.SAFETY_BELT.apply(lambda x: 'Нет' in x)\n",
    "\n",
    "\n",
    "ct = pd.crosstab(data['SAFETY_BELT_NO'], data['target'], normalize=True)\n",
    "ct.plot.bar(stacked=True)\n",
    "plt.legend(title='Погиб')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Странно, но ремень слабо влияет на смертность. Возможно дело в некоректных данных, сложно представить как сотрудники после аварии определяют был ремень или нет. Можно было бы предположить, что только в случае, если человек не выбрался из автомобиля, погиб вероятнее всего. То есть была бы утечка данных таргета. Но нет, распределение не подтверждает, что признак ремень == \"нет\" как-то коррелирует с таргетом. Вероятно, заполняют просто так."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "ct = pd.crosstab(data['osv'], data['target'], normalize=True)\n",
    "\n",
    "ct.plot.bar(stacked=True)\n",
    "plt.legend(title='target')\n",
    "\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Посмотрим на аварийность по регионам с учетом численности населения этих регионов"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# численность населения рф по регионам\n",
    "# http://www.gks.ru/free_doc/new_site/population/demo/Popul2018.xls\n",
    "\n",
    "regions = pd.read_excel('Popul2018.xls', names=['reg', 'value'], header=5, usecols=[0,1])\n",
    "regions['reg'] = regions['reg'].apply(lambda x: x.strip())\n",
    "regions['reg'] = regions['reg']\\\n",
    "    .replace('г. Москва', 'гор. Москва')\\\n",
    "    .replace('г.Санкт-Петербург', 'гор. Санкт-Петербург')\\\n",
    "    .replace('Республика Татарстан', 'Республика Татарстан (Татарстан)')\\\n",
    "    .replace('Республика Адыгея', 'Республика Адыгея (Адыгея)')\\\n",
    "    .replace('Чувашская Республика', 'Чувашская Республика - Чувашия')\\\n",
    "\n",
    "okrug_filter = regions['reg'].str.contains('округ') & ~regions['reg'].str.contains('Чукотский автономный округ')\n",
    "regions = regions[~okrug_filter]\n",
    "\n",
    "d = {elem[0]: elem[1] for elem in regions.values}\n",
    "\n",
    "\n",
    "plt.subplots(figsize=(20,4))\n",
    "ax = sns.barplot(x=regions.reg, y=regions.value)\n",
    "ax.set_xticklabels(ax.get_xticklabels(), rotation=90)\n",
    "ax.set_ylabel(\"Численность населения\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data['region_only'] = data['reg_name'].apply(lambda x: x.split(',')[0])  # выделяем регион ДТП в отдельную колонку\n",
    "\n",
    "by_region = data['region_only'].value_counts()\n",
    "by_region_df = pd.DataFrame(by_region)\n",
    "by_region_df['people'] = list(map(lambda x: d.get(x), by_region.index.values))\n",
    "\n",
    "by_region_df['percent'] = by_region_df['region_only'] / by_region_df['people']\n",
    "\n",
    "\n",
    "plt.subplots(figsize=(20,4))\n",
    "ax = sns.barplot(x=by_region_df.index, y=by_region_df.region_only)\n",
    "ax.set_xticklabels(ax.get_xticklabels(), rotation=90)\n",
    "ax.set_ylabel(\"Число аварий\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.subplots(figsize=(20,4))\n",
    "ax = sns.barplot(x=by_region_df.index, y=by_region_df.percent)\n",
    "ax.set_xticklabels(ax.get_xticklabels(), rotation=90)\n",
    "ax.set_ylabel(\"Отношение числа аварий к населению\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Видно, что в процент аварий к численности населения сильно отличается по регионам, есть более опасные и более безопасные. Возможны пропуски в данных, какой-то из регионов выпал из скрипта (см. Ингушетия), собирающего информацию (особенности кодирования, ошибки в ОКВЕД, предмет отдельного изучения)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Посмотрим на алкоголь"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "data['alco_sum'] = data.ALCO.apply(lambda x: sum(map(int, x)))"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "sns.heatmap(data.corr())"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Есть слабая корреляция алкоголя, погибших и класса дороги"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "zjh7XPrPbsQ7"
   },
   "source": [
    "# 4. Инсайты, найденные зависимости\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Очевидно, сильное значение имеет алкоголь, регион, вид дтп, освещение, тип машины, тип нарушения, в общем почти всё важно. Попробуем загнать все характеристики в модель и посмотрим важность признаков."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "FMUzySiTbsQ8"
   },
   "source": [
    "# 5. Выбор метрики\n",
    "Мы поставили задачу бинарной классификации. Наша задача - получить веса, максимально отделяющие классы друг от друга. Классы не сбалансированы. Есть два основных кандидата на роль метрики при таких обстоятельствах, ROC-AUC и LogLoss. В данном случае нам скорее подходит ROC-AUC, как метрика, максимально отделяющая классы, и потому что нам не так важно насколько именно мы правы или нет (для этого подходит LogLoss). [Чуть более подробные объяснения](https://stats.stackexchange.com/questions/322408/logloss-vs-gini-auc)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "zg00ewoFbsQ9"
   },
   "source": [
    "# 6. Выбор модели\n",
    "\n",
    "Задача бинарной классификации, число признаков, с учетом OneHotEncoding, не превышает 10000 (иначе бы мы были ограничены линейными моделями), поэтому мы можем попробовать несколько алгоритмов и посмотреть, какой из них справится лучше, и уже его настраивать. Мы хотим получить интерпретируемые веса, поэтому не рассматриваем нейронные сети. Вероятно, хорошо подойдут следующие модели:\n",
    "\n",
    "- LogisticRegression\n",
    "- SVC (квадратичная сложность от числа примеров, слишком долго получится на полном датасете, проверим на части)\n",
    "- GradientBoosting\n",
    "- RandomForest\n",
    "- SGDClassifier\n",
    "\n",
    "Конкретные результаты посмотрим после предобработки данных"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "z1BrS1AJbsQ-"
   },
   "source": [
    "# 7. Предобработка данных\n",
    "\n",
    "\n",
    "Мы выбираем для обучения только категориальные признаки, числовых в нашем датасете нет, соответственно нет выбросов, пропуски мы уже заполнили. Раскладываем списки с помощью MultiLabelBinarizer, а простые категории - с помощью LabelBinarizer. Бинарные признаки берем как есть.\n",
    "Отбор признаков осуществим в п.9 плана."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "GLYyT086bsQ_"
   },
   "outputs": [],
   "source": [
    "def prepare_x(multi_list, single_list, binary_list, data):\n",
    "    arr_multicategorial = []\n",
    "    for col in multi_list:\n",
    "        mlb = MultiLabelBinarizer()\n",
    "        transformed = mlb.fit_transform(data[col])\n",
    "        arr_multicategorial.append(pd.DataFrame(transformed, \n",
    "                                                columns=col + '_' + mlb.classes_.astype('object'),\n",
    "                                                ).to_sparse())\n",
    "\n",
    "    arr_categorial = []\n",
    "    for col in single_list:\n",
    "        lb = LabelBinarizer()\n",
    "        transformed = lb.fit_transform(data[col])\n",
    "        arr_categorial.append(pd.DataFrame(transformed,\n",
    "                                           columns=col + '_' + lb.classes_.astype('str').astype('object'),\n",
    "                                           ).to_sparse())\n",
    "\n",
    "\n",
    "    binary_data = data[binary_list]\n",
    "\n",
    "\n",
    "    concated_df = pd.concat(arr_multicategorial+arr_categorial + [binary_data], axis=1).to_sparse()\n",
    "    \n",
    "    return concated_df"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "KCjexihxbsRJ"
   },
   "outputs": [],
   "source": [
    "m = ['t_ts', 'color', 'SAFETY_BELT', 'all_NPDD', 'POL', 't_n', 'OBJ_DTP', 'factor', 'ndu', 's_pog', 'sdor']\n",
    "s = ['DTP_V', 'dor_k', 'dor_z', 'k_ul', 'osv']\n",
    "b = []\n",
    "\n",
    "\n",
    "\n",
    "X = prepare_x(m, s, b, data)\n",
    "y = data['target']  \n",
    "\n",
    "columns_names = X.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "Zhwf68j1bsRO"
   },
   "outputs": [],
   "source": [
    "from sklearn.model_selection import train_test_split, cross_validate\n",
    "\n",
    "# оставляем 10% данных на тест\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.1, random_state=RANDOM_STATE, shuffle=True)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "\n",
    "score_lr = cross_validate(LogisticRegression(), X_train, y_train, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_gbc = cross_validate(GradientBoostingClassifier(), X_train, y_train, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_rf = cross_validate(RandomForestClassifier(), X_train, y_train, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_sgd = cross_validate(SGDClassifier(), X_train, y_train, scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_svc_linear = cross_validate(SVC(kernel='linear', probability=True), X_train[:10000], y_train[:10000], scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_svc_rbf = cross_validate(SVC(probability=True), X_train[:10000], y_train[:10000], scoring='roc_auc')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_lr, score_gbc, score_rf, score_sgd, score_svc_linear, score_svc_rbf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Получается, что несколько алгоритмов показывают хорошие результаты. Можно выделить логистическую регрессию, градиентный бустинг(немного долго) и случайный лес(сильно переобучается, но если над этим посидеть результаты можно улучшить). SVC тоже может дойти до хороших результатов, но квадратичная сложность делает его малоприменимым при нашем датасете\n",
    "\n",
    "Возьмем как основной алгоритм логистическую регрессию."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "kfuAGOQkbsRU"
   },
   "source": [
    "# 8. Кросс-валидация и настройка гиперпараметров модели"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Данных у нас много, не временные ряды, датасет уже произвольно отсортирован при train-test сплите. Поэтому обычный K-Fold."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "ukB5Reh0bsRV"
   },
   "outputs": [],
   "source": [
    "kf = KFold(n_splits=3, random_state=RANDOM_STATE)  # по умолчанию в большинстве методов sklearn"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# обучим модель и посмотрим на веса признаков\n",
    "clf = LogisticRegression()\n",
    "\n",
    "clf.fit(X_train, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "eli5.show_weights(clf, feature_names=columns_names.values, top=50)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Признак SAFETY_BELT при любых значениях, как видно, является утечкой в целевую переменную (и \"Да\", и \"Нет\" имеют высокие положительные веса), потому что его появление может быть объяснимо только в случае, если человек не смог выйти сам из автомобиля, то есть погиб, вылетел из машины или был без сознания.\n",
    "Его, очевидно, лучше исключить.\n",
    "\n",
    "Видно, что мы забыли разложить признак t_n \"технические неисправности\" на составляющие, сделаем это в п.9\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "У логистической регрессии мы можем перебрать следующие параметры в случае бинарной классификации:\n",
    "\n",
    "- penalty : str, ‘l1’ or ‘l2’, default: ‘l2’\n",
    "- C : Сила регуляризации, меньше значения - сильнее регуляризация\n",
    "- class_weight : dict or ‘balanced’, default: None\n",
    "\n",
    "Используем поиск по гриду на части выборки для скорости, С посмотрим приблизительно, более точно в п.10."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "param_grid = dict(penalty=['l1', 'l2'], C=[0.1, 1], class_weight=['balanced', None])\n",
    "gs = GridSearchCV(LogisticRegression(random_state=RANDOM_STATE, n_jobs=-1), param_grid=param_grid, scoring='roc_auc', cv=kf, verbose=10)\n",
    "gs.fit(X_train[:50000], y_train[:50000])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gs.cv_results_"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "plt.title(\"Результаты поиска по гриду\")\n",
    "plt.xlabel(\"Номер набора параметров\")\n",
    "plt.ylabel(\"Score\")\n",
    "plt.ylim(0.78, 0.82)\n",
    "lw = 2\n",
    "plt.plot(range(1, 9), gs.cv_results_['mean_test_score'], label=\"Test score\",\n",
    "             color=\"darkorange\", lw=lw)\n",
    "plt.fill_between(range(1, 9), gs.cv_results_['mean_test_score'] - gs.cv_results_['std_test_score'],\n",
    "                 gs.cv_results_['mean_test_score'] + gs.cv_results_['std_test_score'], alpha=0.2,\n",
    "                 color=\"darkorange\", lw=lw)\n",
    "plt.legend(loc=\"best\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Видно, что гиперпараметры в нашем случае (кроме C, который стоит перебрать) не оказывают существенного влияния, все в пределах стандартного отклонения, но возьмем лучший, и в п.10 настроим C."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gs.best_params_"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "e2W8e5fGbsRa"
   },
   "source": [
    "# 9. Создание новых признаков и описание этого процесса "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "HUCSahjlbsRb"
   },
   "outputs": [],
   "source": [
    "# добавляем временные признаки, дорожные условия, погода, количество машин на дорогах, освещение. Важные признаки\n",
    "data['datetime'] = pd.to_datetime(data['date'] + ' ' + data['Time'], dayfirst=True)\n",
    "\n",
    "data[\"Year\"] = data[\"datetime\"].dt.year\n",
    "data[\"Month\"] = data[\"datetime\"].dt.month\n",
    "data['Weekday'] = data['datetime'].dt.weekday\n",
    "data[\"Hour\"] = data[\"datetime\"].dt.hour\n",
    "data[\"Month_Day\"] = data['datetime'].dt.day\n",
    "data[\"Year_Day\"] = data['datetime'].dt.dayofyear\n",
    "\n",
    "# Year_Day видимо будет переобучением, проверим"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# добавляем нарушения пешехода, упустили в начале\n",
    "def get_sublists_npdd_uch(x):\n",
    "    res = list()\n",
    "    for uch in json.loads(x):\n",
    "        for npdd in uch['NPDD']:\n",
    "            res.append(npdd)\n",
    "    return res\n",
    "\n",
    "data['NPDD_UCH'] = data['uchInfo'].apply(get_sublists_npdd_uch)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# создаем категориальные признаки алкоголя, очевидно одна из основных \n",
    "data['MAX_ALCO'] = data['ALCO'].apply(lambda x: max(map(int, x)) if x else 0).astype(int)\n",
    "data['HIGH_ALCO'] = data['MAX_ALCO'] > 20\n",
    "data['MEDIUM_ALCO'] = (data['MAX_ALCO'] > 10) &  (data['MAX_ALCO'] < 20)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# разворачиваем признак t_n в список, например \"Неисправность рабочей тормозной системы | Не устанавливались\" в две категории\n",
    "def convert_t_n(x):\n",
    "    res = list()\n",
    "    for e in x:\n",
    "        res += e.split('|')\n",
    "    return res\n",
    "        \n",
    "data['t_n'] = data['t_n'].apply(convert_t_n)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "emFDgc2ebsRi"
   },
   "outputs": [],
   "source": [
    "m = ['t_ts', 'color', 'SAFETY_BELT', 'all_NPDD', 'POL', 't_n', 'OBJ_DTP', 'factor', 'ndu', 's_pog', 'sdor', 'NPDD_UCH']\n",
    "s = ['DTP_V', 'dor_k', 'dor_z', 'k_ul', 'osv', 'region_only', 's_pch', 'Month', 'Weekday', 'Hour']\n",
    "b = ['MEDIUM_ALCO', 'HIGH_ALCO']\n",
    "\n",
    "X = prepare_x(m, s, b, data)\n",
    "y = data['target']\n",
    "column_names = X.columns"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# оставляем 10% данных на тест\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.1, random_state=RANDOM_STATE, shuffle=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Оценим результаты с признаком ремня безопасности и без них"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "clf = LogisticRegression(C=0.1,  class_weight='balanced', random_state=RANDOM_STATE)\n",
    "\n",
    "clf.fit(X_train, y_train)\n",
    "score_with_belt = roc_auc_score(y_test,clf.predict_proba(X_test)[:, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def df_without_belt(df):\n",
    "    return df.loc[:, ~df.columns.str.startswith('SAFETY_BELT')]\n",
    "\n",
    "clf = LogisticRegression(C=0.1,  class_weight='balanced', random_state=RANDOM_STATE)\n",
    "clf.fit(df_without_belt(X_train), y_train)\n",
    "score_without_belt = roc_auc_score(y_test,clf.predict_proba(df_without_belt(X_test))[:, 1])"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "score_with_belt, score_without_belt"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Без ремня стало хуже, но, опять таки, по моим ощущениям - это утечка в таргет, невозможно сотруднику в простом ДТП понять, был ремень или нет. Сомнительный признак.\n",
    "Переделываем X без него"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "m = ['t_ts', 'color', 'all_NPDD', 'POL', 't_n', 'OBJ_DTP', 'factor', 'ndu', 's_pog', 'sdor', 'NPDD_UCH']\n",
    "s = ['DTP_V', 'dor_k', 'dor_z', 'k_ul', 'osv', 'region_only', 's_pch', 'Month', 'Weekday', 'Hour']\n",
    "b = ['MEDIUM_ALCO', 'HIGH_ALCO']\n",
    "\n",
    "X = prepare_x(m, s, b, data)\n",
    "y = data['target']\n",
    "column_names = X.columns\n",
    "\n",
    "X_train, X_test, y_train, y_test = train_test_split(X, y, test_size=0.1, random_state=RANDOM_STATE, shuffle=True)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Выбираем 50% самых значимых признаков        "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from sklearn.feature_selection import SelectPercentile\n",
    "\n",
    "sp = SelectPercentile(percentile=50) # выбираем 50% самых значимых признаков\n",
    "\n",
    "X_selected = sp.fit_transform(X_train, y_train)\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_test_selected = sp.transform(X_test)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns = X_test.columns[sp.get_support()]\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "nYRHKhcubsRl"
   },
   "source": [
    "# 10. Построение кривых валидации и обучения"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "GsQbaCyUbsRn"
   },
   "outputs": [],
   "source": [
    "def show_curve(train_scores, valid_scores, param_range, title, xlabel, ymin=0.78, ymax=0.82):\n",
    "    train_scores_mean = np.mean(train_scores, axis=1)\n",
    "    train_scores_std = np.std(train_scores, axis=1)\n",
    "    test_scores_mean = np.mean(valid_scores, axis=1)\n",
    "    test_scores_std = np.std(valid_scores, axis=1)\n",
    "\n",
    "    plt.title(title)\n",
    "    plt.xlabel(xlabel)\n",
    "    plt.ylabel(\"Score\")\n",
    "    plt.ylim(ymin, ymax)\n",
    "    lw = 2\n",
    "    plt.semilogx(param_range, train_scores_mean, label=\"Training score\",\n",
    "                 color=\"darkorange\", lw=lw)\n",
    "    plt.fill_between(param_range, train_scores_mean - train_scores_std,\n",
    "                     train_scores_mean + train_scores_std, alpha=0.2,\n",
    "                     color=\"darkorange\", lw=lw)\n",
    "    plt.semilogx(param_range, test_scores_mean, label=\"Cross-validation score\",\n",
    "                 color=\"navy\", lw=lw)\n",
    "    plt.fill_between(param_range, test_scores_mean - test_scores_std,\n",
    "                     test_scores_mean + test_scores_std, alpha=0.2,\n",
    "                     color=\"navy\", lw=lw)\n",
    "    plt.legend(loc=\"best\")\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "DzEoWaBRbsRq",
    "outputId": "b5195c7f-f66a-403c-f2b0-c5d816e867a4"
   },
   "outputs": [],
   "source": [
    "param_range = np.logspace(-4, 2, 20)\n",
    "train_scores, valid_scores = validation_curve(LogisticRegression(class_weight='balanced', random_state=RANDOM_STATE), X_selected, y_train, \"C\", param_range, n_jobs=-1, cv=kf, scoring='roc_auc', verbose=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "1ELb4Ab9bsRz"
   },
   "outputs": [],
   "source": [
    "show_curve(train_scores, valid_scores, param_range, \"Validation Curve with LogisticRegression\", \"C\")"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "-exiSOzNbsR3"
   },
   "outputs": [],
   "source": [
    "train_sizes, learning_train_scores, learning_valid_scores = learning_curve(LogisticRegression(C=0.1,  class_weight='balanced', random_state=RANDOM_STATE), X_selected, y_train, train_sizes=np.logspace(-2, 0, 10), cv=kf, scoring='roc_auc', verbose=10)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "_BSMFFHQbsR7"
   },
   "outputs": [],
   "source": [
    "show_curve(learning_train_scores, learning_valid_scores, train_sizes, \"Learning Curve with LogisticRegression\", \"Training set size\", 0.75, 0.85)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "train_sizes"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "По кривой валидации видно, что лучшее значение параметра C лежит около 0.1. При дальнейшем увеличении результаты не растут. Кривые обучения близко, ошибка средняя, обусловленная, по всей видимости, случайной природой целевой переменной.\n",
    "\n",
    "Кривая обучения показывает, что результаты на тренировочной и валидационной части сходятся, на полной выборке 500 тысяч строк сойдутся окончательно (На картинке - 2/3 выборки, с учетом KFold3)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "sKxbGYnhbsR_"
   },
   "source": [
    "# 11. Прогноз для тестовой выборки\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# результаты на тестовой выборке\n",
    "clf = LogisticRegression(C=0.1,  class_weight='balanced', random_state=RANDOM_STATE)\n",
    "\n",
    "clf.fit(X_selected, y_train)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "final_score = roc_auc_score(y_test,clf.predict_proba(X_test_selected)[:, 1])\n",
    "final_score"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Результаты на тестовой выборке сравнимы с результатами на кросс-валидации, даже чуть лучше за счет использования всей тренировочной части (на кривых обучения и валидации использовалось 66%)\n",
    "\n",
    "Проведем финальное обучение на полном датасете и посмотрим результаты.\n",
    "\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "X_selected_all = sp.fit_transform(X, y)\n",
    "\n",
    "clf.fit(X_selected_all, y)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "columns_names = X.columns[sp.get_support()]  # получаем имена оставленных колонок\n",
    "\n",
    "eli5.explain_weights(clf, feature_names=columns_names.values, top=200)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab_type": "text",
    "id": "AfePaCnIbsSA"
   },
   "source": [
    "# 12. Выводы"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "dHLl25Z0bsSB"
   },
   "source": [
    "Мы получили интересные результаты, которые могут помочь оценить риск в той или иной ситуации на дороге, выделили признаки на которые стоит обратить внимание при вождении и техническом обслуживании автомобиля.\n",
    "\n",
    "Видно, что выезд на полосу встречного движения является самым страшным с точки зрения последствий нарушением ПДД, пить надо меньше, а пребывание в Краснодарском крае или Воронежской области равновелико по опасности с отсутсвием освещения! Так же стоит проявлять особую бдительность при встрече с крупными транспортными средствами. И аккуратнее ночью!\n",
    "\n",
    "Вдоль дороги лучше вообще не работать или пользоваться всеми возможными бликерами.\n",
    "\n",
    "Неисправность тормозной системы, тонкие шины, разные шины на разных осях - причина гибели!\n",
    "\n",
    "Ну а самое безопасное нарушение - правила пользования общественным транспортом. Ну и на автостоянке сложно погибнуть.\n",
    "\n",
    "Также, в результате работы мы получили интересный качественный датасет, с которым можно работать дальше. К примеру, можно посмотреть зависимость травматичности от марки или модели машины, в случае конкретных участников (мы рассматривалии аварию в целом). Можно посмотреть геоданные, координаты, названия дорог, населенных пунктов для определения наиболее опасных участков.\n"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "colab": {
     "autoexec": {
      "startup": false,
      "wait_interval": 0
     }
    },
    "colab_type": "code",
    "id": "QtQsFHFHbsSG",
    "outputId": "81d8c2c2-7a52-45f4-998d-57239fa98548"
   },
   "source": [
    "Буду рад замечаниям и предложениям\n",
    "\n",
    "\n",
    "[Павел Куликов](http://kulikovpavel.ru)\n",
    "\n",
    "kulikovpavel@gmail.com\n",
    "\n",
    "+7 903 118 37 41"
   ]
  }
 ],
 "metadata": {
  "colab": {
   "default_view": {},
   "name": "project_kulikovpavel.ipynb",
   "provenance": [],
   "version": "0.3.2",
   "views": {}
  },
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.3"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
